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ABSTRACT 

We investigate the spectral properties of the brightest gamma-ray flares of blazars 
detected by the Fermi Large Area Telescope. We search for the presence of spectral 
breaks and measure the spectral curvature on typical time scales of a few days. We 
identify significant spectral breaks in fewer than half of the analyzed flares, but their 
parameters do not show any discernible regularities, and in particular there is no in¬ 
dication for gamma-ray absorption at any fixed source-frame photon energy. More 
interestingly, we find that the studied blazars are characterized by significant spectral 
variability. Gamma-ray flares of short duration are often characterized by strong spec¬ 
tral curvature, with the spectral peak located above 100 MeV. Since these spectral 
variations are observed despite excellent photon statistics, they must reflect tempo¬ 
ral fluctuations in the energy distributions of the emitting particles. We suggest that 
highly regular gamma-ray spectra of blazars integrated over long time scales emerge 
from a superposition of many short-lived irregular components with relatively nar¬ 
row spectra. This would imply that the emitting particles are accelerated in strongly 
turbulent environments. 

Key words: gamma rays: galaxies — quasars: general — quasars: individual: 
3C 454.3, PKS 1510-089 


1 INTRODUCTION 

Blazars, a class of active galactic nuclei, belong to the bright¬ 
est cosmic sources of high-energy (~ GeV) gamma-ray ra¬ 
diation. Their gamma-ray emission is produced by non- 
thermal populations of highly energetic particles in relativis¬ 
tic jets. The mechanism of particle acceleration responsible 
for the emergence of those populations is not yet understood, 
but detailed analysis of the gamma-ray data on blazars has 
the potential to provide more constraints on these processes. 

The Fermi Large Area Telescope (LAT) has unique ca¬ 
pabilities in high-energy astronomy: a very broad spectral 
range (~ 100 MeV — 100 GeV), a very wide field of view 
(~ 60°), and the ability to scan the entire sky every 3 
hours (Atwood et al. 2009). After several years of its mis¬ 
sion, it has collected a vast amount of data on the tem¬ 
poral and spectral behaviour of blazars. Of particular im¬ 
portance are data on the brightest blazars, known as Flat- 
Spectrum Radio Quasars (FSRQs). Their broad-band spec¬ 
tral energy distributions (SEDs) are strongly dominated by 
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the MeV-GeV gamma-ray band (Fossati et al. 1998), and 
they reach gamma-ray fluxes of 10~ 5 ph s -1 cm -2 and higher 
(Abdo et al. 2011). At such high fluxes, they can be analyzed 
at very high temporal and/or spectral resolution. 

The gamma-ray spectra of FSRQs are typically steep, 
with average photon indices of F ~ 2.2 — 2.7 (Abdo et al. 
2010a; Ackermann et al. 2011). However, these spectra are 
not consistent with simple power laws (SPLs); they can be 
better fit by either a log parabola (LP) model, or a broken 
power law (BPL) model. In many blazars, spectral breaks 
are routinely identified, and they can be parametrized by the 
observed photon energy at the break Sbr.obs and the change 
in the photon index AF. In the initial Fermi /LAT data on 
3C 454.3, the brightest blazar of the Fermi era, a break 
was identified in spectra integrated over a month of obser¬ 
vations with .Ebr.obs — 2.4 GeV and Ar ~ 1.2 (Abdo et al. 
2009). In the subsequent studies of 3C 454.3, breaks were 
found with Sbr.oba — 1.0 — 2.8 GeV and Ar ~ 0.6 — 1.0 
(Ackermann et al. 2010; Abdo et al. 2011). In other blazars, 
spectral breaks were identified with break energies ~ 
1.6 — 10 GeV in the source frame (Abdo et al. 2010a). 

It was immediately recognized that such sharp spec- 
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tral breaks cannot be due to a transition to the fast-cooling 
regime of electron energy distribution, which predicts a 
change in photon index of AF = 0.5 (Abdo et al. 2009). For 
a power-law distribution of electrons, it was demonstrated 
that transition to the Klein-Nishina regime of inverse Comp¬ 
ton scattering does not produce a sharp spectral break, 
but rather a smooth cut-off (Tavecchio & Ghisellini 2008; 
Ackermann et al. 2010). However, considering a curved elec¬ 
tron energy distribution, a spectral break can be obtained 
in this regime (Cerruti et al. 2013). The conclusion of 
Abdo et al. (2009) was that these breaks are most likely 
due to a break in the underlying electron energy distribu¬ 
tion. However, Finke & Dermer (2010) noted that this sce¬ 
nario is incompatible with the observed optical/UV spectra 
of 3C 454.3. Instead, they proposed that the observed breaks 
could be explained by a superposition of two spectral compo¬ 
nents, produced by Comptonization of broad emission lines 
and direct radiation of the accretion disk. 

An interesting alternative was proposed by 
Poutanen & Stern (2010), who argued that the gamma-ray 
spectral breaks in blazars could arise from absorption by 
recombination continua of ionized helium (He II). This 
model predicts that the observed breaks should be close 
to f?br — 5 GeV in the source frame. In a subsequent 
study of the gamma-ray spectra of 3C 454.3 at several 
different flux levels, Stern & Poutanen (2011) found a 
weak anti-correlation between the optical depth for He II 
continuum, inferred from AF, and the total gamma-ray 
flux. The implication of strong absorption by ionized helium 
is that gamma-ray emitting regions should be located very 
close (~ 0.1 pc) to the central black hole, possibly shifting 
to larger distances with increasing gamma-ray flux. 

More recent studies have cast doubt on the univer¬ 
sal value of break energies measured in the source frame 
(Harris et al. 2012) and have also suggested that some 
breaks identified in early studies using Fermi /LAT data 
could be artifacts of inaccurate instrument response func¬ 
tions (Harris et al. 2014). Stern & Poutanen (2014) relaxed 
their original claim after updated analysis and now argue 
that only two FSRQs of the nine that they studied display 
spectral breaks at a consistent energy Ebi — 5 GeV, while 
most of these sources show breaks at Ebi — 20 GeV in¬ 
terpreted as due to absorption by hydrogen recombination 
continua. 

Another interesting aspect of the gamma-ray spectra of 
blazars is the spectral curvature, which can be probed by 
fitting log-parabola models to the spectra. Such models re¬ 
turn two parameters — the curvature index /3, and the peak 
photon energy E pea k (see definitions in §2). Typical spec¬ 
tra of FSRQs have Speak < 100 MeV with /? ~ 0.05 — 0.3 
(Ackermann et al. 2011; Harris et al. 2014). It is generally 
thought that the gamma-ray spectra of blazars are sta¬ 
ble, but this problem was not investigated very thoroughly. 
Stern & Poutanen (2011) showed for 3C 454.3 that the cur¬ 
vature index is roughly independent of the gamma-ray flux, 
but that the peak energy is strongly correlated with the 
gamma-ray flux, falling in the range of ~ 10 — 100 MeV. 

In Nalewajko (2013) (hereafter Paper I), a sample was 
selected of the brightest gamma-ray flares of blazars dur¬ 
ing the first four years of the Fermi mission. The study in 
Paper I focused on the temporal properties of the flares, 
such as duration and time asymmetry; however, time varia¬ 


tions of the gamma-ray photon index were also investigated. 
It was demonstrated that many blazar flares exhibit signif¬ 
icant variations of the photon index, and that some very 
short flares (of duration < 1 d) show relatively hard spec¬ 
tra, with T < 2. The most notable example of such flares is 
the MJD 55854 event in PKS 1510-089 (Saito et al. 2013), 
also included in our sample. Additional similar events have 
been identified recently, and they may represent an impor¬ 
tant new class of gamma-ray events observed in luminous 
blazars. 

Here, we analyze the spectral properties of the sample of 
flares selected in Paper I. Rather than attempting to provide 
a general model of the observed gamma-ray spectra, we seek 
to identify and analyze features that are consistent with any 
of the three standard spectral models (SPL, BPL or LP). In 
particular, we focus on two aspects: 1) the occurrence and 
properties of spectral breaks, and 2) spectral curvature and 
spectral variations in the sample. We begin by describing 
our analysis of the spectral properties of the flares in Section 
2. We then present our results on the spectral breaks and 
the spectral curvature in Section 3. This is followed by a 
discussion in Section 4, and conclusions in Section 5. 


2 DATA ANALYSIS 

We investigate the sample of 40 bright blazar flares se¬ 
lected in Paper I (details of the selection procedure are given 
therein). Each flare is a period of time when the observed 
photon flux exceeds half of the peak flux. The minimum 
peak flux for the sample is 7.1 x 10” 6 phs -1 cm -2 . These 
flares have durations 0.5 — 10 days, and they are produced by 
only five blazars: 3C 454.3, PKS 1510-089, PKS 1222+216, 
3C 273, and PKS 0402-362. 

For each flare period, we integrate the binned SED us¬ 
ing the standard analysis tool gtlike from the Fermi /LAT 
ScienceTools software package (v9r27pl). In the analy¬ 
sis, we use the instrument response function P7S0URCE_V6, 
the Galactic diffuse emission model gal_2yearp7v6_v0, the 
isotropic background model iso_p7v6source, events of the 
SOURCE class, the region of interest of radius 10°, and all 
background sources from the 2FGL catalog (Nolan et al. 
2012) within 15°. The energy bins are of equal width and 
uniformly distributed on a logarithmic scale, and they over¬ 
lap with the logarithmic shift equal to 1/3 of the logarithmic 
bin length. In cases of insufficient photon statistics, standard 
2cr upper limits are calculated. The results are shown in Fig 
1 in order of decreasing peak flux (as in Paper I). 

Next, for each flare we performed maximum likelihood 
analysis on the unbinned data in the energy range between 
E m i n = 0.1 GeV and E max = 10 GeV by fitting spectral 
models of an SPL, BPL, and LP. The maximum energy limit 
E max is chosen due to relatively short integration time scales, 
on which only a handful of > 10 GeV photons are detected 
for each flare. 

For the SPL fit, we used the PowerLaw2 spectral model 
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from the Fermi analysis tools, 1 


dJV 

dE 


N{1-T)E~ 
_B 1_r — E 1_ 


(i) 


where dN/dE is the differential photon flux as a function of 
photon energy, T is the spectral index, and -E m in and E max 
are the fixed lower and upper bounds of the energy range. 
Using this model rather than the basic power law fit allowed 
us to avoid artificially selecting a normalization energy and 
instead fix the energy range over which we wished to apply 
the model. In similar fashion, we used the BrokenPowerLaw2 
spectral model for the BPL fit, 


dN ( (E/E hr )~ r ' if E<E hl 

dE \{E/E hl ) 2 otherwise, 


where No is a normalization factor, Ebr is the break energy, 
and Pi and P 2 are the spectral indices in the range where 
E < Ebr and E > Ebr respectively. For the LP fit, we used 
the LogParabola model, 




-(a + ^log(B/B 0 )) 


(3) 


where Eq is a fixed pivot energy, No is the normalization, 
/3 describes the spectral curvature, and a is a parameter 
related to the peak energy, 


Ep ea k 


Eo exp 



(4) 


We set the pivot energy to Eo = 500 MeV, in keeping with 
other analyses (e.g. Harris et al. (2012)), after first testing 
that moving this pivot doesn’t have a significant impact on 
the likelihood value for the fit. 

In analyzing the fit of these three models to the data, 
we did not seek to determine whether the models provided 
a good global description of the observed spectra. Our goal 
was instead to determine which of the three models provided 
the best fit to the data. To do this, we performed an Akaike 
Information Criterion (AIC) test, determining the AIC value 
for each model: 


AIC = 2 k- 21n„Sf, (5) 

where k is the number of free parameters in the model, and 
Jzf is the maximized value of the likelihood function for the 
estimated model. The preferred model is then the one with 
the minimum AIC, and if the difference in the AIC between 
two fits is such that AAIC > 2, the preference for the model 
with the minimum AIC is generally considered to be sta¬ 
tistically significant (refer to Harris et al. (2014), Bozdogan 
(1987), and Lewis et al. (2011) for more information on the 
AIC and AAIC). In our case, if the BPL fit had a AAIC > 2 
relative to the other fits, we declared this spectral break sig¬ 
nificant and recorded the fit parameters in Table 1. 

Uncertainties in most fit parameters were obtained di¬ 
rectly from the fits, with the exception of the parameter 
Ebr • Uncertainties in Ebr were instead estimated from the 
statistical uncertainty corresponding to —2A L = 1 for the 
l-(j confidence region, where L is the log-likelihood function 
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(see Ackermann et al. (2010), Rolke et al. (2005) for discus¬ 
sion of this technique). In cases where Ebr was unbound, we 
rejected the break. 

In many cases, the best-fit BPL model with break en¬ 
ergy Ehri did not yield a significant spectral break, while a 
candidate break could be seen at £#2 # Eb ri- A good exam¬ 
ple of this is flare #10 (3C 454.3, MJD 55294) with i?bri ~ 
0.9 GeV and Ebr 2 ~ 9 GeV. A candidate break at Ebr 2 
may not be identified by the global BPL fit due to the over¬ 
all spectral curvature. This was noted by Stern & Poutanen 
(2011), who found that spectral breaks can be identified 
more robustly when the underlying background spectral 
model is a log-parabola rather than a power-law. In order 
to account for spectral breaks that could have been missed 
by the global ( primary ) BPL fits, for every flare in the sam¬ 
ple we performed an additional search for secondary spectral 
breaks by fitting the BPL, SPL and LP models separately in 
the energy ranges Emin < E < Ebri ( low-energy secondary ) 
and Ebri < E < £ ma x ( high-energy secondary). In the ex¬ 
ample of flare #10, a secondary break was indeed identified 
at Ebr 2 ~ 8.4 GeV. All statistically significant secondary 
spectral breaks are included in Table 1. 


3 RESULTS 

3.1 Spectral breaks 

In analyzing the spectra of the 40 flares, we found that none 
were best fit by a SPL model, 31 were best fit by a BPL 
model, and nine were best fit by a LP model. Of the 31 
BPL-favored spectra, the BPL model was significantly fa¬ 
vored over other models in 15 cases. The breaks from two of 
these cases were rejected because they were unbound, but 
the 13 remaining primary breaks are recorded in Table 1. 
Ten additional statistically significant breaks were found by 
the secondary analysis. They are also recorded in the table 
(labeled to identify them as secondary breaks found either 
on the low- or high-energy side of the primary break of the 
spectrum), for a total of 23 significant breaks detected in the 
40 flares. Figure 1 displays the binned spectrum for each of 
the 40 flares and illustrates the location of any significant 
primary and secondary breaks, as well as the values of the 
spectral index F on either side of the breaks. The primary 
breaks that are not significant are also noted on the spectra 
for reference. It is important to note that the models were 
only fit to the data below 10 GeV; in cases where spectral 
breaks are evident, the observed structure of the spectra 
could extend into the very high energy (VHE) range, but 
this cannot be determined due to very low photon statistics 
above 10 GeV. 

As part of the analysis of the spectral breaks of these 
flares, we sought to test the ‘double absorber’ model put for¬ 
ward by Poutanen & Stern (2010). This model predicts two 
increases in the opacity of the broad-line region to gamma- 
ray-energy photons: one at ~ 5 GeV (in the source frame) 
due to He II recombination, and one at ~ 20 GeV due to H 
I recombination. These opacity increases result in changes 
in the photon index that should be seen clearly as breaks in 
the flare spectra. Because our analysis does not extend be¬ 
yond 10 GeV due to relatively short integration time scales, 
we do not attempt to comment on the presence of the pro¬ 
posed break due to H I, however, we analyze the breaks we 
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found in the 0.1 — 10 GeV range to search for a preference 
for breaks near 5 GeV. In Figure 2, we plot the distribution 
of the redshift-corrected break energies against the peak flux 
of the flare. As can be seen, there is no indication for a break 
preference at 5 GeV or any other fixed energy in the source 
frame. Instead, the break energies seem to be spread uni¬ 
formly between 0.2 GeV < Ebr < 20 GeV. We do not find 
any indication for a correlation between Ebr and other pa¬ 
rameters of the flares, such as peak flux, duration, or time 
asymmetry. 

In Figure 2, we also plot the distribution of the break 
amount, or change in the photon index AF at the break, vs. 
peak flux. For completeness, we also plot the distribution of 
Ar vs. Ebr - For two spectral breaks, we found Ar < 0 (neg¬ 
ative breaks), indicating spectral hardening with increasing 
energy. Both of these breaks are secondary, and they re¬ 
flect local spectral features. For the positive breaks, we find 
that 0.5 < AF < 2. The value of Ar = 0.5 can be under¬ 
stood theoretically as resulting from a break in the electron 
energy distribution N( 7 ) oc i y~ p by A p = 1, which could 
be due to the transition from inefficient to efficient cooling 
(Dermer & Schlickeiser 1993). While a substantial number 
of significant spectral breaks are close to this value in our 
data, roughly the same number of breaks are clearly incon¬ 
sistent with it. The value of AF does not seem to be cor¬ 
related with other parameters of the flare such as observed 
flux, duration or time asymmetry. Stern & Poutanen (2011) 
showed that in the case of 3C 454.3 the value of AT, inter¬ 
preted as a measure of optical depth for absorption from the 
He II continuum, is weakly anticorrelated with the gamma- 
ray luminosity (or flux). We do not find any evidence for 
this in our results. 

3.2 Spectral curvature 

In analyzing the flare spectra, our second goal was to exam¬ 
ine the spectral curvature of all 40 flares, regardless of which 
model provided the best fit to the data. Table 2 contains 
flare details and fit parameters from the LP model applied 
to each of the 40 flares across the full energy range. The spec¬ 
tral curvature parameter f3 spans the range of ~ 0.05 — 0.3, 
which is consistent with the finding that none of the flares 
has a gamma-ray spectrum best described by SPL model. 
It also means that all the gamma-ray spectra are concave, 
which is natural for spectra dominated by a single spectral 
component. Interestingly, we find a very broad range of the 
spectral peak values, with 0.03 GeV < E pea k < 2 GeV. The 
very low values E pea k <0.1 GeV are often insignificant, 
as they are derived from an extrapolation of the LP model 
beyond the observed energy range. Nevertheless, these low 
estimates correspond to the typical shape of the gamma-ray 
spectra of FSRQ blazars, where the unobserved SED peak 
is generally thought to lie somewhere in the ~ 1 — 10 MeV 
range (Fossati et al. 1998; Abdo et al. 2010b). The cases 
where E pea k >0.1 GeV would in general be considered to 
be atypical. 

In Figure 3, we show the distributions of /3 and E pea k 
vs. the flare duration T. We find that flare duration has a 
strong influence on the shape of the gamma-ray spectrum. 
Flares longer than ~ 2.5 d have gently curved spectra with 
f) ~ 0.1 and -Epeak 1 $ 0.1 GeV, whereas shorter flares can 
have a stronger curvature with E pea k > 0.1 GeV. As was 


noted in Paper I, most of the long flares were produced by 
3C 454.3, and they also tend to have more consistent average 
photon indices (r) ~ 2.3, and more symmetric distribution 
of the time asymmetry parameter. On the other hand, the 
short flares are more typical for blazars PKS 1510-089 and 
PKS 1222+216. The latter source stands out by producing 
flares with the highest E pea k values. For completeness, in 
Figure 3 we also show the distribution of /3 vs. E pea k. The 
spectra with E pea k <0.1 GeV by necessity are characterized 
by low values of /3. The highest values of /3 are obtained when 
E pea k ~ 0.2 GeV. 


4 DISCUSSION 

In this work we focus on the gamma-ray spectra of blazars 
integrated during the highest observed gamma-ray fluxes on 
relatively short time scales of T < 10 d. This is a conse¬ 
quence of the flare definition adopted in Paper I, and this ap¬ 
proach distinguishes this work from most studies of gamma- 
ray spectra of blazars that focus on maximizing the photon 
statistics by integrating the spectra on much longer time 
scales (months - years). The results of studies performed 
on longer time scales may not be applicable on shorter time 
scales. A glance at Figure 1 reveals many irregularities that 
are absent in the neat long-term results presented, e.g., by 
Abdo et al. (2010a). In the face of such irregularities, we 
should not expect that these spectra can be well fit by sim¬ 
ple spectral models like SPL, BPL or LP. Instead, we should 
expect at most to produce better or worse approximations 
of the real spectra with these models. 

Of course, one should carefully consider whether these 
irregularities may be due to any statistical or systematical 
errors in the analysis of the Fermi /LAT data. From the 
standard maximum likelihood analysis, we know how many 
reconstructed photons contribute to each flux measurement, 
and so we are not concerned about the photon statistics. The 
systematic errors are less understood, and the Fermi Col¬ 
laboration provides only crude estimates based on the in- 
orbit calibration studies (Ackermann et al. 2012), but these 
values (~ 10 %) are much lower than the observed ampli¬ 
tudes of the spectral fluctuations. That the power density 
spectra of bright blazars are power laws without any breaks 
(Abdo et al. 2010c), including 3C 454.3 in the high state 
(Ackermann et al. 2010), suggests that flux measurements 
are equally accurate at all relevant time scales. We will there¬ 
fore assume that the observed spectral fluctuations are a real 
property of blazars, and not instrumental artifacts. 

Our systematic and unbiased search for the occurrence 
of spectral breaks returned results that can be character¬ 
ized as random. The broad distributions of break energies 
Eb r and break amounts Ar, and the lack of clear corre¬ 
lations with other flare parameters, suggests that there is 
no unique physical mechanism behind them. In the double 
absorber model (Poutanen & Stern 2010), spectral breaks 
should be observed at consistent break energy Ebr — 5 GeV 
in the source frame at all times, unless the gamma-ray ra¬ 
diation is produced far outside the broad-line region. Our 
results do not indicate any preference for this Ebr value, 
which is consistent with the results of Harris et al. (2012). 
While the uncertainty of the intrinsic gamma-ray spectrum 
of blazars may affect the expected position of spectral breaks 
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produced by absorption, it is rather unlikely that this ef¬ 
fect could explain the very broad distribution of Ehi values 
shown in Figure 2. The irregularity of the break parameters 
suggests that they reflect the random spectral fluctuations 
observed in the binned spectra. At much longer integration 
time scales, more regular spectral breaks could arise due to 
non-uniform statistics of such fluctuations. Also, longer in¬ 
tegration times are required in order to address the question 
of ~ 20 GeV breaks arising from absorption by the hydrogen 
Ly-a continuum. 

The main finding of this work is that the spectra of 
long flares (T > 2d) are more regular than the spectra of 
short flares (T < 2 d), which is illustrated by the distri¬ 
bution of the parameters /? and -E pea k of the log-parabola 
fits to the individual spectra (Fig. 3). The short flares of¬ 
ten have their spectral peak within the Fermi /LAT band 
(Speak >0.1 GeV), which is not the case for the long-term 
average spectra of FSRQ blazars (Abdo et al. 2010a). This 
may have profound implications for the theoretical picture of 
dissipation and particle acceleration in relativistic AGN jets. 
The regular (gently broken power-law) gamma-ray spectra 
of blazars observed on long time scales may generally be 
superpositions of many simple components. Each of those 
components may have a narrowly peaked particle energy 
distribution and produce a short subflare contributing to 
the overall observed light curve. If the emission of such a 
component is coherent in photon energy and time, it may 
also be coherent in space. In short, we propose that the emit¬ 
ting regions of blazars are macroscopically turbulent, with 
individual eddies producing short radiation spikes of narrow 
spectrum. 

Of course, it is not a new idea that blazar jets are fun¬ 
damentally turbulent (Begelman et al. 1984; Jones 1988). 
Observations that the power density spectra of blazar light 
curves are consistent with simple power laws (Abdo et al. 
2010c) provide strong evidence of this. Dedicated studies 
have been done in order to explain the observed blazar light 
curves in terms of emission produced in relativistic tur¬ 
bulence (Marscher 2014; Calafut & Wiita 2014). What we 
demonstrate here is that Fermi /LAT is able to resolve these 
turbulent fluctuations in energy when observations are inte¬ 
grated on time scales T < 2 d. If the maximum size of the 
turbulent eddies is limited by the jet radius R at distance 
r, the observed time scale of a single fluctuation is limited 
by t v ar < R/{T>c) ~ r/(F 2 c) ~ 3 d (r/1 pc)(r/20) -2 , where 
we assumed that R = dr, 6 ~ 1/F is the jet opening angle, 
F ~ F is the Doppler factor, and F is the jet Lorentz factor. 
If the gamma-ray radiation in blazars is produced mainly at 
distance scale r ~ 0.1 — 1 pc (Nalewajko et al. 2014), this 
time scale would be broadly consistent with our results. 

Macroscopic turbulence may be generated in relativis¬ 
tic strongly-magnetized jets by current-driven instability 
(CDI). Global 3D general relativistic simulations of black- 
hole jets suggest that they can be marginally stable to 
global CDI modes (McKinney & Blandford 2009). For suf¬ 
ficient dissipation efficiency, the most desirable situation is 
to have a jet stable globally and unstable locally (Begelman 
1998). Recent numerical results demonstrate that this is the¬ 
oretically feasible (Porth & Komissarov 2014). Other mech¬ 
anisms that can produce macroscopic turbulence in rel¬ 
ativistic jets include plasma heating in the jet bound¬ 
ary layer resulting from interaction with the environment 


(Kohler & Begelman 2015), or magnetic reconnection gen¬ 
erating relativistic bulk motions and/or relativistic heating 
(Giannios et al. 2009). 

If the broad spectra of blazars are a superposition of 
many narrow components that can be resolved on daily time 
scales, there is no need to explain the regular long-term av¬ 
erage spectra by microscopic processes like diffusive shock 
acceleration at relativistic shock waves (Spitkovsky 2008). 
Instead, one can consider a process that energizes many par¬ 
ticles to a similar energy, like relativistic magnetic reconnec¬ 
tion. Recent numerical studies of relativistic reconnection 
(Sironi & Spitkovsky 2014; Guo et al. 2014; Werner et al. 
2014) demonstrate that for sufficiently high magnetization 
a > 10 , the particles form a very hard power-law distri¬ 
bution N( 7 ) oc 7 _p with p < 2, therefore, most energy is 
contained (and can be radiated) at the high-energy end of 
the distribution defined by the total energy content in the 
system. 

The gamma-ray spectra of flaring blazars presented 
here can be used to constrain the underlying energy dis¬ 
tribution of ultra-relativistic electrons, under the assump¬ 
tion that the gamma-ray radiation is produced by inverse 
Compton scattering of external photon fields (Sikora et al. 
2009). However, the relation between the energy distribu¬ 
tion of electrons and the spectrum of radiation scattered 
by them is complicated; for example, even in the Thomson 
regime the curvature parameter for the electrons is about 
four times greater than the curvature parameter for the ra¬ 
diation (Massaro et al. 2006). Knowledge of the simultane¬ 
ous synchrotron radiation spectra in the NIR/optical/UV 
band and detailed modeling of the broad-band blazar spec¬ 
tra would be necessary in order to determine the electron 
energy distribution. Our proposition for the turbulent ori¬ 
gin of the gamma-ray spectra of blazars further complicates 
this problem by requiring a multi-zone approach. 

In Paper I, a dichotomy was revealed between the tem¬ 
poral properties of the brightest gamma-ray flares of blazars. 
On one hand, most flares produced by 3C 454.3 are long, 
with complex light curves (multiple subflares of compara¬ 
ble peak flux), without clear time asymmetry. On the other 
hand, most flares produced by PKS 1510-089 and PKS 
1222+216 are short, with simple light curves, and a ten¬ 
dency for the flux decay time scale to be longer than the 
flux raising time scale. This dichotomy was suggested to be 
an observational effect, with the viewing angle of the jet 
much smaller in the case of 3C 454.3, resulting in the more 
uniform Doppler beaming of all emitting regions. Now we 
add to this picture the systematic differences between the 
observed gamma-ray spectra. The spectra of long flares are 
more regular than the spectra of short flares, as the for¬ 
mer consist of more elementary narrow components. The 
interpretation of the source dichotomy in terms of the view¬ 
ing angle is consistent with this, as more uniform Doppler 
beaming (smaller viewing angle) is required to observe more 
spectral components at comparable flux levels (hence more 
regular integrated spectra). 

Two of the sources contributing to the sam¬ 
ple of flares studied here, PKS 1222+216 and PKS 
1510-089, were detected in the VHE gamma-ray band 
(above 70 GeV) by ground-based Cherenkov tele¬ 
scopes. The MAGIC observation of PKS 1222+216 at 
MJD 55364.9 (Aleksic et al. 2011) coincided with the be- 
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ginning of Flare # 12 (Tanaka et al. 2011). The H.E.S.S. 
observations of PKS 1510-089 at MJD 54910-54918 
(H.E.S.S. Collaboration et al. 2013) coincided with Flare 
#20, and MAGIC observations of PKS 1510-089 between 
MJD 55975 and MJD 56020 (Aleksic et al. 2014) coincided 
with a series of flares: Flares #19, #29, #28, and #32. 
The observed vF v flux of PKS 1222+216 at 100 GeV was 
~ 10 -10 ergs -1 cm -2 , while that of PKS 1510-089 was 

much lower, ~ 10 -11 ergs -1 cm -2 . The Fermi/ LAT spec¬ 
tra of the flares coinciding with the VHE detections do not 
appear to be different from the spectra of other flares for 
the same sources. Unfortunately, because of the very limited 
photon statistics above 10 GeV at short integration times, 
we cannot determine which of these spectra extend to the 
VHE band at detectable levels. However, by extrapolating 
the best-fit log-parabola models (and ignoring the problem 
of possible absorption of the VHE gamma-rays by intrinsic 
and external soft radiation fields), we can make an educated 
guess which of the flares were most likely to be detected at 
100 GeV at the high level of vF v > 10 -10 ergs -1 cm -2 : 
for PKS 1222+216, Flares #8, #12, #26 and #35; for 
PKS 1510-089, Flare #11; and for 3C 454.3, Flare #1. 
PKS 1222+216 appears to be most likely to be detected in 
the VHE band due to its consistently hard gamma-ray spec¬ 
trum (Nalewajko 2013), although so far — to our knowledge 
— it was detected only once. In the case of 3C 454.3, at red- 
shift 2 = 0.859, any VHE detection would be extraordinar¬ 
ily significant in constraining the extragalactic background 
light (EBL), and we suggest that it should be vigorously 
pursued by the current (H.E.S.S., MAGIC and VERITAS) 
and future ( CTA ) VHE observatories. 


5 CONCLUSIONS 

We performed a spectral analysis with the Fermi /LAT of 
the sample (selected in Paper I) of the 40 brightest gamma- 
ray flares of blazars (FSRQs) detected in the first four years 
of the Fermi mission. The gamma-ray spectra are integrated 
over relatively short time scales T < 10 d, and they show 
significant and variable departures from the long-term av¬ 
erage spectra of the same sources. We performed a uniform 
search for the occurrence of spectral breaks. The break en¬ 
ergies show a broad distribution and no preference for the 
fixed value of 5 GeV in the source frame predicted by the 
double-absorber model of Poutanen & Stern (2010). In or¬ 
der to compare the basic structures of the observed spec¬ 
tra, we fitted them with a log parabola model and found 
an interesting trend of the model parameters with the flare 
duration. Short flares (T < 2 d) often show a strong spec¬ 
tral curvature with the SED peak within the Fermi/ LAT 
range E pea ^ > 0.1 GeV, while all long flares show a mild 
spectral curvature with a SED peak below the LAT range. 
The dichotomy between typical properties of flares observed 
in sources 3C 454.3 vs. PKS 1510-089 and others, first de¬ 
scribed in Paper I, is extended to include differences between 
the observed gamma-ray spectra. The extrapolation of the 
log-parabola spectral models to 100 GeV identifies several 
candidates for VHE detections. 

We suggest that the irregular gamma-ray spectra ob¬ 
served by the Fermi/ LAT for short blazar flares reflect 
macroscopic turbulence in relativistic jets with individual 


emitting regions (turbulent eddies) having a narrow en¬ 
ergy distribution of emitting particles. The superposition of 
many such spectral components peaking at different energies 
would then result in the regular power-law spectra observed 
over long time scales. This would eliminate the need for lo¬ 
calized stochastic particle acceleration in the blazar zones 
of relativistic AGN jets. Further investigation is necessary 
in order to understand the underlying energy distribution of 
ultra-relativistic electrons. 
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Figure 1. Time-integrated flare spectra displaying the i/F u flux (ergs -1 cm -2 ) vs observed energy (GeV) of the 40 brightest gamma-ray 
blazar flares detected by Fermi/ LAT. Statistically significant primary spectral breaks (those identified with a BPL fit to the entire energy 
range 0.1 — 10 GeV) and secondary spectral breaks (those identified with a BPL fit to only the energy range above or below the primary 
break energy) are indicated by broken lines that illustrate the mean value of the spectral index T on either side of the break. Primary 
breaks that are not statistically significant are indicated by open circles. Dotted lines show the best-fit primary log-parabola models. 
The color of the spectrum indicates the host blazar. 
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Figure 2. Distributions of the source-frame break energy E\^ r 
(GeV), the change in photon index Ar = T 2 — Ti, and the 
peak photon flux F pea k (10 —5 phs -1 cm -2 ), plotted for all iden¬ 
tified spectral breaks in the top 40 Fermi /LAT gamma-ray flares. 
The shape and color of each point indicate the host blazar. Solid 
symbols indicate statistically significant primary and secondary 
breaks, whereas unfilled symbols indicate primary breaks that 
were not statistically significant. 


Figure 3. Distributions of the spectral curvature (3 and the spec¬ 
tral peak position F pea ^ (resulting from a log parabola fit to the 
spectrum) versus the duration T (days) of the flare, plotted for 
the 40 brightest Fermi /LAT gamma-ray flares. The shape and 
color of each point indicate the host blazar. 
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Table 1. Statistically significant breaks discovered in the spectra of the 40 brightest Fermi /LAT gamma-ray flares. Breaks were recorded 
if the spectrum was better-fit by a BPL rather than SPL or LP model, with a AAIC > 2 for the BPL compared to the next-best fit. 
MJDp ea k is the moment of flux peak, -F pe ak is the peak photon flux in units of 10 —6 phs -1 cm -2 , T is the duration of the flare in days, 
£»br (obs) and E \> r (source) are the break energies, in GeV, in the observer and source frame respectively, and F i and F 2 are the spectral 
indices on either side of the break energy. The break classification identifies the break as resulting from a fit in the entire energy range 
0.1 — 10 GeV (“primary”) or resulting from a fit only in the energy range below or above the primary break (“secondary: low” and 
“secondary: high”). 


# 

Blazar 

MJD peak 

Break classification 


Fpeak 

T 

E hl (obs) 

Fbr (source) 



r 1 



r 2 

i 

3C 454.3 

55520.0 

secondary: 

low 

77.2 

± 

2,4 

3.9 

0 35+ 0 05 
U.OO—006 

0.65l°:?? 

2.03 

± 

0.03 

2.25 

± 

0.04 

1 

3C 454.3 

55520.0 

secondary: 

high 

77.2 

± 

2,4 

3.9 

6 46 +0 ‘ 62 
o.^to _2 07 

12.00±|;“ 

2.43 

± 

0.10 

3.59 

± 

0.31 

4 

3C 454.3 

55550.3 

secondary: 

high 

23.6 

± 

1.3 

8.3 

0 - 19 + 3.79 
O.lZ—2 47 

15.09ll;°| 

2.54 

± 

0.05 

3.32 

± 

0.43 

5 

3C 454.3 

55167.8 

secondary: 

low 

21.8 

± 

1.5 

4.3 

0 15+ 001 
U, 1 O - 0.02 

C4 ^ 

0 q 
d d 
+ 1 
00 

CN 

cd 

1.62 

± 

0.27 

2.22 

± 

0.05 

6 

3C 454.3 

55567.8 

secondary: 

high 

19.2 

± 

1.6 

5.4 

c 00+1-81 

o.04_2 22 

10 . 81 I 413 

2.36 

± 

0.05 

3.06 

± 

0.38 

7 

3C 454.3 

55301.5 

primary 


16.4 

± 

1.8 

3.5 

1 64 +0 ' 51 
1,o ^-0.32 

3.051!°! 

2.25 

± 

0.04 

2.89 

± 

0.20 

8 

PKS 1222+216 

55316.6 

secondary: 

high 

15.6 

± 

0.9 

0.8 

7 16 +1,58 

1 ,lo -0.80 

10.26ll°! 

1.97 

± 

0.17 

3.46 

± 

0.64 

10 

3C 454.3 

55294.1 

secondary: 

high 

15.3 

± 

0.8 

9.5 

8 35+ 112 
0 . 00-2 68 

15.521°°® 

2.59 

± 

0.09 

4.13 

± 

0.81 

12 

PKS 1222+216 

55365.8 

primary 


14.2 

± 

1.0 

2.3 

U.OO-005 

0.80lJ°® 

1.76 

± 

0.08 

2.27 

± 

0.06 

12 

PKS 1222+216 

55365.8 

secondary: 

low 

14.2 

± 

1.0 

2.3 

0 17+ 001 
u,± - 0.02 

0 24 -1 " 0 01 
0.03 

2.76 

± 

0.36 

1.58 

± 

0.13 

13 

3C 454.3 

55163.1 

primary 


11.8 

± 

1.4 

3.0 

1 O7-+0.27 
-0.28 

9 cc+0-51 
Z.OO-052 

2.17 

± 

0.06 

3.19 

± 

0.23 

14 

3C 454.3 

55305.5 

primary 


11.1 

± 

1.5 

1.9 

u.oo_ 0>06 

1 03+ 0 - 41 
i.UO_o 10 

1.98 

± 

0.09 

2.94 

± 

0.15 

20 

PKS 1510-089 

54917.0 

primary 


9.9 

± 

0.7 

2.4 

1 [ 50”*" 0 ' 30 
i.OU_o 28 

2 04+ 0 - 41 
z,u ^_0.38 

2.03 

± 

0.05 

2.75 

± 

0.17 

21 

PKS 1510-089 

55851.9 

primary 


9.9 

± 

1.0 

0.7 

1 QK + 0- 25 

i.uo_o.29 

1 42 - *" 0,33 
x -^ z -0.39 

2.04 

± 

0.15 

3.24 

± 

0.48 

22 

3C 454.3 

55195.2 

primary 


9.7 

± 

0.8 

2.2 

i-971°;« 

o.OO-o 56 

2.12 

± 

0.06 

3.68 

± 

0.46 

25 

3C 454.3 

55327.2 

primary 


8.8 

± 

1.2 

3.1 

1-80^:33 

3 34+ 0 - 66 
°-°^-0.62 

2.23 

± 

0.06 

2.98 

± 

0.26 

26 

PKS 1222+216 

55342.1 

secondary: 

high 

8.7 

± 

0.8 

1.6 

1 93+ 146 

J..yo_o 71 

9 iyiy-\-2.09 
• ' — 1.02 

1.82 

± 

0.08 

2.24 

± 

0.14 

30 

PKS 1510-089 

54961.8 

primary 


8.2 

± 

0.7 

0.6 

0-S6l°;“ 

0 7fi+ 0 - 30 
u -' 0.15 

2.04 

± 

0.24 

3.56 

± 

0.49 

33 

3C 454.3 

55154.8 

primary 


7.8 

± 

0.9 

2.4 

O 

to 

'I 

1 + 
O O 

b b 

0 c;n+ 0 08 

U.OU-o 07 

1.79 

± 

0.20 

2.50 

± 

0.09 

35 

PKS 1222+216 

55377.5 

primary 


7.6 

± 

0.7 

1.5 

U.OO — o 23 

o-5ilg;g 

1.86 

± 

0.17 

2.33 

± 

0.10 

37 

PKS 1510-089 

55876.1 

secondary: 

high 

7.4 

± 

0.8 

1.7 

0 51+ 0 - 53 
u,ol -0.10 

o.83l°;I° 

4.54 

± 

0.90 

2.31 

± 

0.16 

38 

PKS 1222+216 

55234.0 

primary 


7,4 

± 

0.8 

1.1 

0 21 + 003 
u,z± —0.02 

0 29”*” 0 ' 05 
u.zy_o 02 

0.50 

± 

0.59 

2.49 

± 

0.13 

40 

3C 454.3 

55091.6 

primary 


7.1 

± 

0.8 

5.5 

3 09+ 0,76 
o.u»_0 84 

3 75 - *" 1 - 41 

2.35 

± 

0.05 

4.63 

± 

0.97 
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Table 2. Log-parabola fit parameters for the brightest gamma-ray flares of blazars. MJD pea k is the moment of flux peak, -Fpeak is the 
peak photon flux in units of 10 —6 phs -1 cm -2 , t\ is the flux doubling time scale, £2 is the flux halving time scale, T = ( t\ +£ 2 ) is the flare 
duration, a is the spectral index at the pivot energy, /3 is the spectral curvature, and -Speak is the implied energy of the observer-frame 
vF v spectral peak in units of GeV. All times are in units of days. 


# 

Blazar 

MJD peak 


-^peak 

tl 

£2 

T 



a 



P 


-^peak 

1 

3C 454.3 

55520.0 

77.2 

± 

2.4 

2.9 

1.0 

3.9 

2.20 

± 

0.01 

0.09 

± 

0.01 

0.17 

± 0.02 

2 

3C 454.3 

55526.9 

33.3 

± 

1.8 

1.1 

3.9 

5.0 

2.31 

± 

0.02 

0.13 

± 

0.02 

0.15 

± 0.03 

3 

PKS 1510-089 

55853.8 

26.6 

± 

2.0 

0.2 

0.4 

0.6 

1.76 

± 

0.08 

0.23 

± 

0.05 

0.85 

± 0.18 

4 

3C 454.3 

55550.3 

23.6 

± 

1.3 

3.2 

5.1 

8.3 

2.27 

± 

0.02 

0.11 

± 

0.01 

0.14 

± 0.02 

5 

3C 454.3 

55167.8 

21.8 

± 

1.5 

1.2 

3.1 

4.3 

2.29 

± 

0.03 

0.11 

± 

0.02 

0.14 

± 0.04 

6 

3C 454.3 

55567.8 

19.2 

± 

1.6 

2.9 

2.5 

5.4 

2.25 

± 

0.02 

0.08 

± 

0.02 

0.09 

± 0.04 

7 

3C 454.3 

55301.5 

16.4 

± 

1.8 

1.3 

2.2 

3.5 

2.32 

± 

0.04 

0.07 

± 

0.03 


<0.10 

8 

PKS 1222+216 

55316.6 

15.6 

± 

0.9 

0.4 

0.4 

0.8 

1.69 

± 

0.05 

0.13 

± 

0.03 

1.67 

± 0.56 

9 

PKS 1510-089 

55872.8 

15.3 

± 

1.2 

0.2 

0.7 

0.9 

2.00 

± 

0.06 

0.14 

± 

0.04 

0.51 

± 0.11 

10 

3C 454.3 

55294.1 

15.3 

± 

0.8 

5.9 

3.6 

9.5 

2.36 

± 

0.02 

0.08 

± 

0.01 

0.05 

± 0.02 

11 

PKS 1510-089 

55867.8 

14.2 

± 

1.5 

0.5 

1.5 

2.0 

2.11 

± 

0.06 

0.09 

± 

0.04 

0.27 

± 0.12 

12 

PKS 1222+216 

55365.8 

14.2 

± 

1.0 

1.0 

1.3 

2.3 

1.97 

± 

0.04 

0.09 

± 

0.03 

0.60 

± 0.14 

13 

3C 454.3 

55163.1 

11.8 

± 

1.4 

2.0 

1.0 

3.0 

2.30 

± 

0.05 

0.13 

± 

0.04 

0.16 

± 0.07 

14 

3C 454.3 

55305.5 

11.1 

± 

1.5 

0.8 

1.1 

1.9 

2.36 

± 

0.06 

0.22 

± 

0.05 

0.22 

± 0.05 

15 

PKS 1510-089 

55746.2 

11.1 

± 

1.4 

0.4 

0.5 

0.9 

2.26 

± 

0.09 

0.08 

± 

0.07 


<0.25 

16 

PKS 1510-089 

54948.0 

10.6 

± 

2.3 

0.9 

1.3 

2.2 

2.35 

± 

0.06 

0.15 

± 

0.05 

0.16 

± 0.07 

17 

3C 273 

55095.3 

10.3 

± 

0.9 

0.4 

1.6 

2.0 

2.39 

± 

0.06 

0.07 

± 

0.05 


<0.09 

18 

3C 273 

55090.5 

10.2 

± 

0.9 

0.4 

1.2 

1.6 

2.36 

± 

0.07 

0.21 

± 

0.06 

0.21 

± 0.06 

19 

PKS 1510-089 

55980.6 

10.0 

± 

1.3 

0.5 

0.8 

1.3 

2.11 

± 

0.08 

0.28 

± 

0.07 

0.41 

± 0.06 

20 

PKS 1510-089 

54917.0 

9.9 

± 

0.7 

1.9 

0.5 

2.4 

2.12 

± 

0.04 

0.10 

± 

0.03 

0.28 

± 0.08 

21 

PKS 1510-089 

55851.9 

9.9 

± 

1.0 

0.2 

0.5 

0.7 

2.26 

± 

0.12 

0.16 

± 

0.10 

0.22 

± 0.14 

22 

3C 454.3 

55195.2 

9.7 

± 

0.8 

0.6 

1.6 

2.2 

2.23 

± 

0.06 

0.16 

± 

0.04 

0.24 

± 0.06 

23 

PKS 1222+216 

55310.7 

9.6 

± 

0.8 

0.5 

0.7 

1.2 

1.73 

± 

0.08 

0.21 

± 

0.05 

0.96 

± 0.24 

24 

3C 454.3 

55323.5 

9.1 

± 

0.9 

4.2 

1.8 

6.0 

2.36 

± 

0.04 

0.08 

± 

0.03 

0.05 

± 0.04 

25 

3C 454.3 

55327.2 

8.8 

± 

1.2 

1.9 

1.2 

3.1 

2.29 

± 

0.05 

0.08 

± 

0.04 

0.08 

± 0.08 

26 

PKS 1222+216 

55342.1 

8.7 

± 

0.8 

0.8 

0.8 

1.6 

1.78 

± 

0.05 

0.08 

± 

0.03 

1.93 

± 1.16 

27 

3C 273 

55202.9 

8.7 

± 

1.1 

0.3 

0.6 

0.9 

2.54 

± 

0.12 

0.33 

± 

0.11 

0.22 

± 0.07 

28 

PKS 1510-089 

55990.8 

8.5 

± 

0.7 

5.2 

1.1 

6.3 

2.26 

± 

0.03 

0.08 

± 

0.02 

0.10 

± 0.05 

29 

PKS 1510-089 

55983.1 

8.4 

± 

0.8 

0.5 

0.4 

0.9 

2.10 

± 

0.09 

0.24 

± 

0.07 

0.41 

± 0.08 

30 

PKS 1510-089 

54961.8 

8.2 

± 

0.7 

0.3 

0.3 

0.6 

2.62 

± 

0.18 

0.32 

± 

0.17 

0.19 

± 0.11 

31 

3C 454.3 

55214.3 

8.0 

± 

0.9 

0.8 

0.8 

1.6 

2.37 

± 

0.09 

0.23 

± 

0.08 

0.22 

± 0.08 

32 

PKS 1510-089 

56002.4 

7.8 

± 

1.0 

1.8 

0.9 

2.7 

2.40 

± 

0.06 

0.10 

± 

0.04 

0.07 

± 0.06 

33 

3C 454.3 

55154.8 

7.8 

± 

0.9 

0.5 

1.9 

2.4 

2.30 

± 

0.06 

0.09 

± 

0.05 

0.10 

± 0.10 

34 

PKS 1510-089 

55767.6 

7.6 

± 

0.8 

0.4 

0.5 

0.9 

1.98 

± 

0.09 

0.17 

± 

0.07 

0.54 

± 0.14 

35 

PKS 1222+216 

55377.5 

7.6 

± 

0.7 

0.4 

1.1 

1.5 

2.15 

± 

0.07 

0.05 

± 

0.04 


<0.26 

36 

3C 454.3 

55283.9 

7.6 

± 

1.1 

1.3 

1.4 

2.7 

2.30 

± 

0.07 

0.15 

± 

0.05 

0.19 

± 0.08 

37 

PKS 1510-089 

55876.1 

7.4 

± 

0.8 

0.7 

1.0 

1.7 

2.28 

± 

0.07 

0.05 

± 

0.05 


<0.12 

38 

PKS 1222+216 

55234.0 

7.4 

± 

0.8 

0.7 

0.4 

1.1 

2.20 

± 

0.11 

0.18 

± 

0.09 

0.29 

± 0.12 

39 

PKS 0402-362 

55827.5 

7.3 

± 

1.0 

0.3 

1.8 

2.1 

2.18 

± 

0.08 

0.19 

± 

0.07 

0.32 

± 0.09 

40 

3C 454.3 

55091.6 

7.1 

± 

0.8 

4.4 

1.1 

5.5 

2.41 

± 

0.04 

0.07 

± 

0.03 


<0.07 
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